In [4]:
percorsoIn = "/home/protoss/Documenti/TESI/thesis/codici/matlabbo/dopplero.mat"
struttura = scipy.io.loadmat(percorsoIn)['infos']
freqCorrOrig = struttura['freqniu'][0,0][0]
velPerPosOrig = struttura['velperpos'][0,0][0]
print(freqCorrOrig, velPerPosOrig)
In [1]:
import tensorflow as tf
import numpy
import scipy.io
import time
#nb le fft da 8192 sono intervallo freq da x Hz a y Hz
# da 4096 sono intervallo freq da x Hz a y Hz
percorsoIn = "/home/protoss/Documenti/TESI/thesis/codici/matlabbo/quadHWI52.mat"
#percorsoIn = "/home/protoss/Documenti/TESI/thesis/codici/matlabbo/quadrato4096Quad.mat"
quadrato = scipy.io.loadmat(percorsoIn)['quad']
#print(quadrato)
#carico file dati
percorsoDati = "/home/protoss/Documenti/TESI/DATI/dati9mesi52HWI.mat"
#percorsoDati = "/home/protoss/Documenti/TESI/DATI/datiunmese.mat"
struttura = scipy.io.loadmat(percorsoDati)['job_pack_0']
tempi = struttura['peaks'][0,0][0]
frequenze = struttura['peaks'][0,0][1]
pesi = (struttura['peaks'][0,0][4]+1)
#nb: picchi ha 0-tempi
# 1-frequenze
# 2-pesi
#headers vari
securbelt = 4000 #WARNING DA DEFINIRE MEGLIO SECURBELT
tFft = 8192
stepFreq = 1/tFft
enhancement = 10
stepFreqRaffinato = stepFreq/enhancement
nstepsFreq = securbelt+(numpy.amax(frequenze)-numpy.amin(frequenze) + stepFreq + 2*stepFreqRaffinato)/stepFreqRaffinato
#print(nstepsFreq)
#tempi
#da definire epoca come mediana di tempi di tutto il run
epoca = struttura['basic_info'][0,0]['run'][0,0]['epoch'][0,0][0,0]
#per doppler corr
veloc = struttura['basic_info'][0,0]['velpos'][0,0][0:3,:]
nTempi = struttura['basic_info'][0,0]['ntim'][0,0][0,0]
primoTempo = struttura['basic_info'][0,0]['tim0'][0,0][0,0]
indices = struttura['basic_info'][0,0]['index'][0,0][0]
#print(indices.shape, veloc.shape, nTempi)
In [105]:
nTempi = (numpy.amax(tempi)-numpy.amin(tempi))*3600*24/4096
nFrequenze = (numpy.amax(frequenze)-numpy.amin(frequenze))/stepFreq
print(nTempi, nFrequenze)
In [76]:
#"vanilla"
quadratoNP = quadrato[0]
start = time.time()
#nPunti = quadrato.shape[0]
velocitas = numpy.zeros((3,indices[nTempi]-1))
print(velocitas.shape)
for i in numpy.arange(0,nTempi-1):
intervallo = numpy.arange(indices[i]-1,indices[i+1]-1)
for j in intervallo:
velocitas[:,j] = veloc[:,i]
velPerPosIndex = numpy.dot(quadratoNP,velocitas)
#print(velPerPos)
frequenzeIndex = frequenze / (1+velPerPosIndex)
time.time()-start
#for punto in numpy.arange(nPunti):
# velPerPos = numpy.dot(quadrato[punto],velocitas)
# frequenzeIndex = frequenze / (1+velPerPos)
Out[76]:
In [9]:
#un for tolto
quadratoNP = quadrato[0]
start = time.time()
#nPunti = quadrato.shape[0]
velocitas = numpy.zeros((3,frequenze.size))
for i in numpy.arange(0,nTempi-1):
inizio = indices[i]-1
fine = indices[i+1]-1
velocitas[:,inizio:fine] = veloc[:,i:i+1]
velPerPosIndex = numpy.dot(quadratoNP,velocitas)
divisoreIndex = 1+velPerPosIndex
#print(velPerPos)
frequenzeUnfor = frequenze / divisoreIndex
time.time()-start
#for punto in numpy.arange(nPunti):
# velPerPos = numpy.dot(quadrato[punto],velocitas)
# frequenzeIndex = frequenze / (1+velPerPos)
Out[9]:
In [116]:
#un for tolto
quadratoNP = quadrato[0]
indicesOpt = indices-1
inizi = indicesOpt[:-1]
fini = indicesOpt[1:]
def primaFor():
velocitas = numpy.zeros((3,frequenze.size))
for i in numpy.arange(0,nTempi-1):
velocitas[:,inizi[i]:fini[i]+1] = veloc[:,i:i+1]
velPerPosIndex = numpy.dot(quadratoNP,velocitas)
divisoreIndex = 1+velPerPosIndex
freq = frequenze / divisoreIndex
return freq
def primaProd():
velPerPosIndex = numpy.dot(quadratoNP,veloc)
velocitas = numpy.zeros(frequenze.size)
for i in numpy.arange(0,nTempi-1):
velocitas[inizi[i]:fini[i]+1] = velPerPosIndex[i:i+1]
divisoreIndex = 1+velocitas
freq = frequenze / divisoreIndex
return freq
%time freqPrimaFor = primaFor()
%time freqPrimaProd = primaProd()
#for punto in numpy.arange(nPunti):
# velPerPos = numpy.dot(quadrato[punto],velocitas)
# frequenzeIndex = frequenze / (1+velPerPos)
In [92]:
print(frequenze.size)
quadratoNP = quadrato[0]
starts = indicesOpt[:-1]
ultimiStarts = numpy.where(starts ==frequenze.size)[0]
starts = starts[1:ultimiStarts[0]]
ends = indicesOpt[1:]
ultimiEnds = numpy.where(ends ==frequenze.size)[0]
ends = ends[1:ultimiEnds[0]+1]-1
velPerPosIndex = numpy.dot(quadratoNP,veloc)
def slices2arr_gather(ar, starts, ends, N):
id_arr = numpy.zeros((N),dtype=int)
id_arr[starts[1:]] = 1
c = id_arr.cumsum()
id_arr[ends[1:]] =id_arr[ends[1:]]-1
out = numpy.where(id_arr.cumsum()==0, 0, ar[c])
out[starts[0]:ends[0]] = ar[0]
out[ends] = out[ends-1]
return out
%time velocitas = slices2arr_gather(velPerPosIndex, starts, ends, frequenze.size)
divisoreIndex = 1+velocitas
frequenzeVect = frequenze / divisoreIndex
In [16]:
matrix = numpy.random.randint(2, size = 100).astype(float).reshape(10,10)
times = numpy.nonzero(matrix)[0]+1
positions = numpy.nonzero(matrix)[1]
indexes = numpy.where(numpy.diff(times)>0)[0]+1
indexes = numpy.concatenate(([0],indexes, [times.size]))
nTimes = numpy.size(indexes)-1
speeds = numpy.random.rand(nTimes)
def slices2arr_for():
expandedSpeeds = numpy.zeros(positions.size)
for i in numpy.arange(0,nTimes):
expandedSpeeds[starts[i]:ends[i]] = speeds[i]
return expandedSpeeds
def slices2arr_gather(ar, starts, ends, N):
id_arr = numpy.zeros((N),dtype=numpy.int64)
id_arr[starts[1:]] = 1
c = id_arr.cumsum()
%time id_arr[ends[1:]] =id_arr[ends[1:]]-1
out = numpy.where(id_arr.cumsum()==0, 0, ar[c])
out[starts[0]:ends[0]] = ar[0]
out[ends] = out[ends-1]
return out
def slices2arr_addat(ar, starts, ends, N):
id_arr = numpy.zeros((N),dtype=numpy.int64)
id_arr[starts[1:]] = 1
c = id_arr.cumsum()
%time numpy.add.at(id_arr, ends[1:],-1)
out = numpy.where(id_arr.cumsum()==0, 0, ar[c])
out[starts[0]:ends[0]] = ar[0]
out[ends] = out[ends-1]
return out
starts = indexes[:-1]
ends = indexes[1:]
%time a= slices2arr_for()
starts = indexes[:-1]
ends = indexes[1:]-1
%time b= slices2arr_gather(speeds, starts, ends, positions.size)
%time c=slices2arr_addat(speeds, starts, ends, positions.size)
In [19]:
print(numpy.nonzero(a-b))
In [249]:
Out[249]:
In [6]:
sessione = tf.Session()
tempiTF = tf.constant(tempi, dtype=tf.float64)
velTF = tf.constant(veloc, dtype=tf.float64)
freqTF = tf.constant(frequenze, dtype=tf.float64)
quadTF = tf.constant(quadrato, dtype=tf.float64)
#tempiTF = tempiTF -tf.reduce_min(tempiTF)
tempiTF = tempiTF - primoTempo
tempiUnici, indiciUniciRipetuti, conteggi = tf.unique_with_counts(tempiTF)
indiciUniciRipetutiTF = sessione.run(indiciUniciRipetuti)
indiciUnici = tf.cumsum(tf.pad(tf.unique_with_counts(indiciUniciRipetuti)[2],[[1,0]]))[:-1]
nIndici = tf.size(tempiUnici)
differenze = tempiUnici[1:]-tempiUnici[:-1]
stepTempi = tf.reduce_min(differenze)
# il primo elemento nella peakmap è a tempi successivi del primo elemento delle velocità
tempiPrima = tf.reduce_min(tempiUnici)/stepTempi
tempiPrima = tf.cast(tempiPrima, dtype=tf.int64)
# vedo dove sono i salti
stepSaltati = differenze/stepTempi
stepSaltati = tf.cast(stepSaltati, dtype=tf.int64)
# ora seleziono solo gli indici dell'array di velocità
# ai quali il relativo step temporale ha frequenze non nulle nella peakmap
indiciVel = tf.zeros([tf.size(tempiUnici)-1],dtype=tf.int64) + tempiPrima
indiciVel = indiciVel + tf.cumsum(stepSaltati)
indiciVel = tf.concat([[tempiPrima], indiciVel],0)
# qui va messa parte per calcolare i pesi, che per ora non mettiamo
#def doppcorr(i):
#qui faccio la doppl corr
velPerPos = tf.matmul(tf.reshape(quadTF[0], [1,3]),velTF)
velPerPos = velPerPos
velPerPos = velPerPos[0]
velPerPosTF = sessione.run(velPerPos)
#nel programma originale c'è un prodotto elementwise infattibile
#(tra due array di diversa dimensione)
# c'è anche un for che considera buchi di tempi dentro e prima
#qui cerco tutti i salti temporali vedendo quando tra un elemento e l'altro
#c'è una differenza maggiore (un multiplo) dello step minimo
# e quindi definisco un nuovo array di velocità con i soli i valori
# che mi interessano
#velUtili = velPerPos[indiciVel]
##WARNING GATHER SBALLA IL VALORE DI POCHISSIMO MA ABBASTANZA PER ROVINARE CORR DOPPLER
velUtili = tf.gather(velPerPos, indiciVel)
velUtiliTF = sessione.run(velUtili)
divisore = tf.gather_nd(velUtili,indiciUniciRipetuti[:,None])
divisoreTF = sessione.run(divisore)
freqCorr = tf.divide(freqTF,(divisore +1))
%time freqCorrTF = sessione.run(freqCorr)
In [31]:
print(veloc.shape)
print(numpy.unique(tempi).size)
tempiTot = numpy.arange(primoTempo, numpy.amax(tempi),stepTempi)
print(tempiTot.size)
Out[31]:
In [126]:
print(veloc.shape)
print(numpy.arange(primoTempo*24*3600, numpy.amax(tempi)*24*3600, 4096).size)
print(numpy.arange(numpy.amin(tempi)*24*3600, numpy.amax(tempi)*24*3600, 4096).size)
In [118]:
start = time.time()
tempiDC = numpy.round(((tempi-numpy.amin(tempi))*3600*24)/4096).astype(numpy.int64)
#tempiDC = ((tempi-primoTempo)*3600*24)#.astype(numpy.int64)
tempiUnici,indiciUnici,indiciOrig = numpy.unique(tempiDC, return_index=True,return_inverse=True)
print(tempiUnici)
# qui va messa parte per calcolare i pesi che per ora non mettiamo
#qui faccio la doppl corr
nIndici = numpy.size(tempiUnici)
#differenze = numpy.round(numpy.diff(tempiUnici)).astype(numpy.int64)
differenze = numpy.diff(tempiUnici)
stepTempi = numpy.amin(differenze)
# il primo elemento nella peakmap è a tempi successivi del primo elemento delle velocità
#tempiPrima =(numpy.amin(tempiUnici)/stepTempi).astype(numpy.int64)
# vedo dove sono i salti
stepSaltati = numpy.round(differenze/stepTempi).astype(numpy.int64)
#stepSaltati = differenze/stepTempi
#print(stepSaltati)
#stepSaltati = numpy.concatenate(([tempiPrima],stepSaltati))
#stepSaltati = numpy.concatenate(([0],stepSaltati))
indiciVel = numpy.zeros(numpy.size(tempiUnici))+1 #+tempiPrima
indiciVel[1:] = indiciVel[1:]+numpy.cumsum(stepSaltati)
indiciVel = numpy.round(indiciVel).astype(numpy.int64)
# e quindi definisco un nuovo array di velocità con i soli i valori
# che mi interessano
velPerPos = numpy.dot(quadrato,veloc)
velPerPos = velPerPos
velPerPos = velPerPos[0]
velUtili = velPerPos[indiciVel]
#velUtili = numpy.take(velPerPos, indiciVel)
# infine faccio la divisione
divisore = velUtili[indiciOrig]+1
#divisore = numpy.take(velUtili, indiciOrig)+1
freqCorrNP = frequenze/divisore
stop = time.time()
print(stop-start)
from matplotlib import pyplot
%matplotlib notebook
pyplot.figure(figsize=(10, 8))
a = pyplot.scatter(numpy.arange(freqPrimaFor.size),freqPrimaFor-freqCorrNP, s =10)
In [94]:
from matplotlib import pyplot
%matplotlib notebook
pyplot.figure(figsize=(10, 8))
a = pyplot.scatter(numpy.arange(velPerPos.size),velPerPosIndex-velPerPos, s = 10)
In [22]:
from matplotlib import pyplot
%matplotlib notebook
pyplot.figure(figsize=(10, 8))
a = pyplot.scatter(numpy.arange(tempi.size),tempi, s = 10)
In [93]:
import numpy
times = numpy.array([0.1, 0.1, 0.2, 0.2])
positNP = numpy.array([58.98962402, 58.9921875, 60.00390625, 60.00878906])
vel = numpy.array([0.99705114,0.99974157])
uniqueTimes = numpy.unique(times, return_index=True)
uniqueIndices = uniqueTimes[1]
uniqueTimes = uniqueTimes[0]
numIndices = numpy.size(uniqueTimes)-1
iterator = numpy.arange(numIndices)+1
for i in iterator:
positNP[uniqueIndices[i-1]:uniqueIndices[i]] = posit[uniqueIndices[i-1]:uniqueIndices[i]]*vel[i-1]
In [97]:
import tensorflow as tf
times = tf.constant([0.1, 0.1, 0.2, 0.2])
positTF = tf.constant([58.98962402, 58.9921875, 60.00390625, 60.00878906])
vel = tf.constant([0.99705114,0.99974157])
uniqueTimes, uniqueIndicesRepeated, counts = tf.unique_with_counts(times)
uniqueIndices = tf.cumsum(tf.pad(tf.unique_with_counts(uniqueIndicesRepeated)[2],[[1,0]]))[:-1]
out = tf.gather_nd(vel,uniqueIndicesRepeated[:,None])*positTF
sessione = tf.Session()
outpos = sessione.run(out)
In [65]:
import numpy
times = numpy.array([0.1, 0.1, 1.5, 1.9, 1.9, 1.9])
posit = numpy.array([2.1, 3.5, 0.4, 1.3, 2.7, 3.5])
vel = numpy.array([0.5,0.7,1.0])
uniqueTimes = numpy.unique(times, return_index=True)
uniqueIndices = uniqueTimes[1]
uniqueTimes = uniqueTimes[0]
numIndices = numpy.size(uniqueTimes)-1
iterator = numpy.arange(numIndices)+1
for i in iterator:
posit[uniqueIndices[i-1]:uniqueIndices[i]] = posit[uniqueIndices[i-1]:uniqueIndices[i]]*vel[i-1]
In [98]:
print(positNP)
print(outposTF)
print(positNP-outposTF)
In [152]:
nFrequenze = (1/stepFrequenza).astype(int)
print(nFrequenze)
tempi = tempi-primoTempo
# qui va messa parte per calcolare i pesi, che per ora non mettiamo
#qui faccio la doppl corr
velPerPos = numpy.dot(quadrato,veloc)
velPerPos = velPerPos + 1
velPerPos = velPerPos[0]
#nel programma originale c'è un prodotto elementwise infattibile
#(tra due array di diversa dimensione)
# c'è anche un for che lavora sull'array index nella struttura basicinfo
# ma array index non considera buchi di tempi dentro né prima, solo dopo
#frequenze = tf.multiply(frequenze, velPerPos)
print(numpy.size(velPerPos))
In [154]:
tempiUnici = numpy.unique(tempi, return_index=True)
indiciUnici = tempiUnici[1]
tempiUnici = tempiUnici[0]
In [155]:
nIndici = numpy.size(tempiUnici)-1
indici = numpy.arange(nIndici)+1
differenze = numpy.zeros(nIndici)
for i in indici:
differenze[i-1] = tempiUnici[i]-tempiUnici[i-1]
stepTempi = numpy.amin(differenze)
print(stepTempi)
print(stepTempi*nTempi)
print(numpy.amax(tempi)/stepTempi)
print(numpy.size(tempiUnici))
#print(sessione.run(frequenze))
tempiVeri = numpy.arange(nTempi)*stepTempi
tempiPrima = (numpy.amin(tempi)/stepTempi).astype(int)
tempiInmezzo = ((numpy.amax(tempi)-numpy.amin(tempi))/stepTempi).astype(int)
tempiDopo = nTempi-tempiInmezzo-tempiPrima
print(tempiPrima, tempiDopo, tempiInmezzo)
In [151]:
stepSaltati = (differenze/stepTempi).astype(int)
indiciVel = numpy.zeros(numpy.size(tempiUnici))
indiciVel[0] = tempiPrima
indiciVel[1:numpy.size(indiciVel)] = tempiPrima+numpy.cumsum(stepSaltati)
indiciVel = indiciVel.astype(numpy.int32)
#for i in numpy.arange(numpy.size(stepSaltati)):
# indiciVel[i] = tempiPrima+1+i
velUtili = velPerPos[indiciVel]
for i in indici:
frequenze[indiciUnici[i-1]:indiciUnici[i]] = frequenze[indiciUnici[i-1]:indiciUnici[i]]/velUtili[i]
print(numpy.size(frequenze))
In [ ]:
In [61]:
#c'è qualche problema con le divisioni con numpy
import tensorflow as tf
times = tf.constant([0.1, 0.1, 1.5, 1.9, 1.9, 1.9])
posit = tf.constant([2.1, 3.5, 0.4, 1.3, 2.7, 3.5])
vel = tf.constant([0.5,0.7,1.0])
uniqueTimes, uniqueIndicesRepeated, counts = tf.unique_with_counts(times)
uniqueIndices = tf.cumsum(tf.pad(tf.unique_with_counts(uniqueIndicesRepeated)[2],[[1,0]]))[:-1]
#out = tf.gather_nd(vel,uniqueIndicesRepeated[:,None])*posit
sessione = tf.Session()
#outpos = sessione.run(out)
print(sessione.run(tf.gather_nd(vel,uniqueIndicesRepeated[:,None])))